#include "structCellBase.hpp"
#include "meshTypes.hpp"
#include "structMeshFunctions.hpp"
#include "exceptions.hpp"
#include <cassert>

namespace hms {

template<typename CellType>
StructCellBase<CellType>::StructCellBase( const MeshType& mesh, Index index ) :
	StructCellBase( mesh, index, toDoubleIndex(mesh, index) )
{}

template<typename CellType>
StructCellBase<CellType>::StructCellBase(
	const MeshType& mesh, Index x, Index y
) :
	StructCellBase( mesh, toSingleIndex(mesh, x, y), {x,y} )
{}

template<typename CellType>
StructCellBase<CellType>::StructCellBase(
	const MeshType& mesh, Index index, std::pair<Index, Index> xy
) :
	m_mesh  {mesh},
	m_index {index},
	m_i     {xy.first},
	m_j     {xy.second}
{}


template<typename CellType>
auto StructCellBase<CellType>::i() const -> Index { return m_i; }

template<typename CellType>
auto StructCellBase<CellType>::j() const -> Index { return m_j; }

template<typename CellType>
auto StructCellBase<CellType>::meshImp() const -> const MeshType& {
	return m_mesh;
}

template<typename CellType>
auto StructCellBase<CellType>::indexImp() const -> Index { return m_index; }

template<typename CellType>
auto StructCellBase<CellType>::edgesImp() const -> EdgeIndices {
	assert( !this->isGhost() &&
		"StructCellBase::edges(Index) is for domain cells only!"
		"For ghost cells, use Cells::edge(Index) instead!" );

	Index nHorizontalEdges { m_mesh.nCellsX()*(m_mesh.nCellsY()+1) };
	return {
		m_i + m_j*m_mesh.nCellsX(),
		m_i + m_j*(m_mesh.nCellsX()+1) + nHorizontalEdges + 1,
		m_i + (m_j+1)*m_mesh.nCellsX(),
		m_i + m_j*(m_mesh.nCellsX()+1) + nHorizontalEdges
	};
}

template<typename CellType>
auto StructCellBase<CellType>::edgeImp() const -> Index {
	assert( this->isGhost() &&
		"StructCellBase::edge(Index) is for ghost cells only!"
		"For domain cells, use Cells::edges(Index) instead!");

	Index nHorizontalEdges { m_mesh.nCellsX()*(m_mesh.nCellsY()+1) };
	return
		m_j < 0 ? /* bottom ghost */
			m_i :
			m_i < 0 ? /* left ghost */
				nHorizontalEdges + m_j*(m_mesh.nCellsX()+1) :
				m_i >= m_mesh.nCellsX() ? /* right ghost */
					nHorizontalEdges + (m_j+1)*(m_mesh.nCellsX()+1) - 1 :
					m_i + nHorizontalEdges - m_mesh.nCellsX(); /* top ghost */
}

template<typename CellType>
auto StructCellBase<CellType>::nEdgesImp() const -> char {
	if ( !this->isGhost() ) return 4;
	else return 1;
}

template<typename CellType>
auto StructCellBase<CellType>::isBoundaryImp() const -> bool {
	/* this check is still required, because the check below is only valid
	 * for domain cells */
	if ( this->isGhost() ) return false;
	return
		m_i == 0 || m_i == m_mesh.nCellsX() - 1 ||
		m_j == 0 || m_j == m_mesh.nCellsY() - 1;
}

template<typename CellType>
auto StructCellBase<CellType>::domainCellImp() const -> CellType {
	assert( this->isGhost() && "domainCell can only be called on ghost cells!" );
	Index i, j;
	if ( m_i < 0 ){
		i=0;
		j=m_j;
	} else if ( m_j < 0 ){
		i=m_i;
		j=0;
	} else if ( m_i >= m_mesh.nCellsX() ){
		i=m_mesh.nCellsX()-1;
		j=m_j;
	} else {
		i=m_i;
		j=m_mesh.nCellsY()-1;
	}
	return m_mesh.cell(i,j);
}


template<typename CellType>
auto StructCellBase<CellType>::isGhostImp() const -> bool {
	return
		m_i < 0 ||
		m_j < 0 ||
		m_i >= m_mesh.nCellsX() ||
		m_j >= m_mesh.nCellsY();
}

template<typename CellType>
auto StructCellBase<CellType>::neighboursImp() const -> ArrayXi {
	/* with simplified indexing, only ghost cells need special treatment */
	if ( !this->isGhost() ){
		static constexpr Index ng {1};
		return Array4i{
			m_index - (m_mesh.nCellsX() + 2*ng),
			m_index + 1,
			m_index + (m_mesh.nCellsX() + 2*ng),
			m_index - 1
		};
	}
	/* bottom ghost cells */
	if ( m_j < 0 )
		return Array1i{
			toSingleIndex(m_mesh, m_i, 0)
		};
	/* left ghost cells */
	if ( m_i < 0 )
		return Array1i{
			toSingleIndex(m_mesh, 0, m_j)
		};
	/* right ghost cells */
	if ( m_i >= m_mesh.nCellsX() )
		return Array1i{
			toSingleIndex(m_mesh, m_mesh.nCellsX()-1, m_j)
		};
	/* top ghost cells */
	return Array1i{
		toSingleIndex(m_mesh, m_i, m_mesh.nCellsY()-1)
	};
}

template<typename CellType>
auto StructCellBase<CellType>::verticesImp(bool) const -> ArrayXi {
	/* bottom left index for domain cells */
	Index i { m_j*(m_mesh.nCellsX()+1) + m_i };
	if ( !this->isGhost() ) {
		/* corners 0, 1, 2, 3 (counter-clockwise) */
		return Array4i{ i, i+1, i+m_mesh.nCellsX()+2, i+m_mesh.nCellsX()+1 };
	} else {
		return m_mesh.edge( this->edge() ).vertices();
	}
}

template<typename MeshType>
auto toSingleIndex( const StructMeshBase<MeshType>& mesh, Index i, Index j )
-> Index
{
	assertValidCellIndex(mesh,i,j);
	static constexpr Index ng {1};
	return i + ng + (j + ng) * (mesh.nCellsX() + 2*ng);
}

template<typename MeshType>
auto toDoubleIndex( const StructMeshBase<MeshType>& mesh, Index cell )
-> std::pair<Index, Index>
{
	assertValidCellIndex(mesh,cell);
	static constexpr Index ng {1};
	return {
		cell % ( mesh.nCellsX()+2*ng ) - ng,
		cell / ( mesh.nCellsX()+2*ng ) - ng
	};
}


template class StructCellBase<UniCell>;
template class StructCellBase<RectCell>;
template class StructCellBase<StructCell>;

#define INSTANTIATE_FOR_STRUCTURED_MESHTYPES(MeshType) \
template Index toSingleIndex( const StructMeshBase<MeshType>&, Index x, Index y ); \
template std::pair<Index, Index> toDoubleIndex( const StructMeshBase<MeshType>&, Index );

#include "structMeshInstantiator.hpp"

}